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Abstract 

We discuss the influence of periodic orbits on the dissociation of a model diatomic molecule 
driven by a strong bichromatic laser fields. Through the stability of periodic orbits we analyze 
the dissociation probability when parameters like the two amplitudes and the phase lag between 
the laser fields, are varied. We find that qualitative features of dissociation can be reproduced 
by considering a small set of short periodic orbits. The good agreement with direct simulations 
demonstrates the importance of bifurcations of short periodic orbits in the dissociation dynamics 
of diatomic molecules. 
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I. INTRODUCTION 



The dissociation behavior of molecules driven by bichromatic fields with commensurate 
frequencies has emerged as a rich research subject, especially for the control of molecular 
processes |2[ . The interplay of the two radiation fields opens up many new dissociation 
pathways, and it is well known that the relative phase between the two fields can affect 
these pathways drastically, so much so that the relative phase can be used as a means to 
control the outcome of the reaction 3, |4J . However, the mechanisms by which the relative 
phase controls the dissociation behavior are less well-known. The relative phase is a very 
convenient control parameter since it does not require additional energy input from the fields 
(as opposed to their amplitudes). 

The two-color laser-driven dissociation of molecules is also of great interest to researchers 
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ias been extensive in the past three decades 
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20j. Among these is Ref. 



3, 



7J where 



the dissociation probability of a model diatomic molecule exposed to a two-color laser field 
was investigated for various parameters using direct simulations of the classical-mechanical 
equations. 

In this paper, we report how the dissociation probability obtained in Ref. 7J by direct 
numerical simulations can be understood qualitatively using a linear stability analysis of a 
small set of periodic orbits. Our main result is that for most parameter values the principal 
features of the dissociation probability can be reproduced using two short periodic orbits 
(with the period equal to the one of the fields), and in particular by the identification of 
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the main bifurcations which have a drastic effect on the dissociation probability. In this 
way, our approach allows the qualitative prediction, with significant time savings, of the 
dynamics as parameters are varied. Indeed, the typical time necessary for the computation 
of a periodic orbit and its stability is of the order of its period which is also the period of 
the field. Our findings echo similar ones obtained in the microwave ionization of Rydberg 
atoms in a strong bichromatic field for which qualitative agreement has been obtained with 
experimental data and a quantitative agreement with quantal simulations 21]. 



After describing the model in Sec. [TTJ we briefly summarize the method which monitors 



221 ] in Sec. IHII In Sec. HV] we analyze the 



the position and the residue of periodic orbits 
dynamics using a selection of short periodic orbits, as parameters are varied in Sec. II V Al 
We relate a linear stability measure (the residue of a given periodic orbit) to the dissociation 
probability in Sec. IIVBI Good agreement is found for most parameter values. A discrepancy 
is observed for small phases and is discussed in Sec. IIVCI 

II. THE MODEL 

The Hamiltonian (in atomic units) of a diatomic molecule exposed to a strong bichromatic 
field with a phase lag can be modeled as 

2 

H(r, p, r) = — V D\l — e~ a ^)] 2 

2m L J 

+ (r - r e ) [A x sin(ftir) + A 2 sin(fi 2 r + </>)], (1) 

where the parameters are the reduced mass m, the dissociation energy D, and the equilibrium 
distance r e . Here we assume that the envelopes of the pulses are constant since the pulse 
duration has a very minor impact on this system, as suggested in Ref. 7|. The relevant 
dimensionless variables are f = a(r — r e ), p = p /V2Dm, t = a^/2D/mr, F { = Ai/(2Da), 



uji — Oj/ \j2Da 2 /m. In these new coordinates, the Hamiltonian ([T]) is 

H(r,P,t) = P l+ l -(l-e-") 2 

+r {F x smuxt + F 2 sm{u 2 t + 4>)) , (2) 

where f and p are canonically conjugate. 

In what follows, we model the hydrogen fluoride (HF) molecule in which m = 1732, 
D = 0.2101, r e = 1.75 and a = 1.22 (all in atomic units.) We consider a field with two 
commensurate frequencies such that uoi = u 2 /3 = 0.28. We notice that Hamiltonian ([2]) is 
time-periodic with period 2ir/ui. 



III. RESIDUE METHOD 



The general idea of the residue method is to follow a set of periodic orbits as parameters 
are varied in order to determine qualitative properties of the dynamics. As it was shown 



in other, simi 
dynamics 21 



ar problems, short periodic orbits play the role of organizing centers for the 
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24|. Higher-order periodic orbits give more refined details of the dynamics, 
especially on longer time scales. In principle, for atomic and molecular systems where 
short pulses are considered, only short periodic orbits should influence the dynamics. We 
determine the location of a periodic orbit (given by its number of intersections with the apt 
section) using a modified Newton-Raphson multi-shooting algorithm as 



Poincare surface o: 

described in Ref. 25j. The initial conditions from which the Newton map are iterated is 
determined in two possible ways : By a quick inspection of the Poincare section (which is 
easier if the periodic orbit is elliptic since it shows resonant islands around the periodic orbit 
considered) or by continuation of periodic orbits for other values of parameters (which is 
clearly the optimal way since periodic orbits usually deform continuously as parameters are 



varied, and even bifurcate). We also monitor the linear stability properties of these periodic 
orbits which are obtained by integrating the reduced tangent flow along the periodic orbit 



Tr _,2£V/~ ~ ,\ rrt 



dt 



where 



/ \ 

1 



\- 10 / 



and V H is the two-dimensional Hessian matrix (composed of second 



derivatives of H with respect to its canonical variables f and p). The initial condition 
is J~° = I 2 (the two-dimensional identity matrix). For a periodic orbit with period T 
where T = 2im/uji (n being the number of intersections with the Poincare section), the 
two eigenvalues of the monodromy matrix J T ( which make the pair (A, 1/A)) determine 
the stability properties. The determinant of J T is equal to 1 since the flow is volume 
preserving. If the spectrum is (e luJ , e~ luJ ) , the periodic orbit is elliptic (stable, except in some 
rare cases); or hyperbolic if the spectrum is (A, 1/A) with A e K* (unstable). Through the 



use of Greene's residue R 
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the stability properties can be deduced in a concise form. If R e]0, 1[, the periodic orbit is 
elliptic; if R < or R > 1 it is hyperbolic; and if R = and R = 1, it is parabolic. 

For a given periodic orbit, we follow its location in phase space and its residue as pa- 
rameters are varied. There are three parameters in our problem: Two amplitudes F\ and 
F 2 , and a phase lag (p. We compute R(Fi, F 2 , (ft) and identify the points in parameter space 
where bifurcations occur. We are interested in the bifurcations whenever a periodic orbit 
is likely to change its linear stability, which occurs in particular for R(Fi, F 2 , (ft) = or 
R(F±, F 2 , 0) = 1. In general such bifurcations (based on a linear stability analysis) will also 
play an important role in the nearby phase space region (by continuity in phase space). 
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IV. DISSOCIATION PROBABILITY 



A. Identification of fundamental periodic orbits 

Figure [T] shows a Poincare section (stroboscopic plot of phase space with period 2-k/uji) 
of Hamiltonian ([2]) for amplitudes F\ = 0.18, F 2 = 0.02 and phase lag = 0. We notice 
that an elliptic island is present at the entrance of the dissociation channel. At the center of 
this island sits an elliptic periodic orbit with one intersection with the Poincare surface of 
section (i.e. with period 2tt/ui). Standard Hamiltonian dynamics show that the trajectories 
that are likely to dissociate can become trapped around the resonant island for a while 
before finding an escape route. Therefore, this particular periodic orbit, which we call O e , 
plays a crucial role in the dissociation probability, and is the focus of the current paper : We 
investigate its role as the parameters (Fx, F 2 , 0) are varied. Due to the symmetry (0 1— > 7r — 0) 
the fundamental domain of variations of is [0,7r[. We also restrict the amplitudes to 



We anticipate two factors which can influence dissociation: One is the location of this 
orbit, and the other is a change of its stability. In Fig. [2J we represented the position of O e 
or, more precisely, its action and angle variables as defined by 10] 



where E(r,p) = p 2 + (1 — e _r ) 2 , for a typical set of parameters Fi = 0.18 and F 2 = 0.02, 
while is varied. The main conclusion is that the variation of the action and angle of O e do 
not seem to be linked to the variation of the dissociation probability since they do not vary 
significantly as is increased. In addition they are not monotonic functions of in contrast 



(F 1 ,F 2 ) e [0,0.22] x [0,0.06]. 
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with the dissociation probability. Hence, the position of the specific periodic orbit does not 
appear to play a significant role (at least in the range of parameters considered). 

The second possible mechanism based on a bifurcation has a more drastic influence on 
dissociation. In order to monitor this bifurcation properly, we also need to follow the asso- 

n 

dated hyperbolic orbit [22]. To see an increase of dissociation we need to ensure that the 
hyperbolic orbit stays hyperbolic while the elliptic one turns hyperbolic (in order to discard 
a stability exchange which would not affect the dissociation probability significantly). A 
typical residue plot (the residue as a function of 0) is shown in Fig. [3] for F\ = 0.18 and 
F 2 = 0.02. We notice that an increase of is always associated with an increase of the 
residue. At = 0.75, the residue crosses unity and a bifurcation (which is a period dou- 
bling) occurs. This increase of hyperbolicity is associated with increased chaos and hence 
more dissociation. This is in agreement with the direct simulations of Ref. 7j. 

We should point out that the computation of Ref. 7J uses initial conditions in the ground 
state (E = 0.045), which is lower than the set of periodic orbits we consider (E = 0.29) in 
Fig.[H This justifies the importance of the chosen periodic orbits for dissociating trajectories. 
For other values of parameters we consider, the energy level of the periodic orbit O e is always 
well above E = 0.045. 

B. Residue contour plots in parameter space 

For fixed values of 0, we vary the amplitudes F\ and F 2 and compute the residue values of 
the periodic orbit O e . We depict the contour plots of these residues in Fig. @]m the (F\, F 2 ) 
plane for = 0, n/6, it/ 2 and it. These figures should be compared to Fig. 2 of Ref. [7J. In 
Fig. HI darker regions represent lower residue values, and hence are expected to reveal more 
stable dynamics ( smaller dissociation probability) for the corresponding parameters. 



On these residue contours we have superposed the curves R(Fi, F 2 ) = 1 (white dashed 
curves) through which the periodic orbit O e changes its stability. The corresponding bifur- 
cation from elliptic to hyperbolic linear stability indicates a possible increase in dissociation 
probability due to an increase of hyperbolicity in this region of phase space. Since the 
trajectories are no longer trapped when O e is hyperbolic, it is expected that the dissoci- 
ation probability increases. This is in good agreement with Ref. [7| where we notice that 
for 0=7r/6, 7r/2 and 7r, there is a qualitative agreement in the shape of the dissociation 
probability. In particular, the following features are reproduced for a fixed value of : 
the non-monotonicity for = 7r/6 as F 2 is increased (with a fixed value of F±), and the 
monotonicity for = 7r/2 and 7r with a sharper downward bifurcation curve (dashed line) 
for = tt as F 2 is increased in the region F2 G [0,0.05]. Our analysis also confirms that 
the stabilization effect decreases when is increased (from to 71). For = case, the 
contour plot shows agreement with direct simulations for most regions in the (F 1; F 2 ) plane 
once more, and the property that it reproduces the two upper-right bumps observed in the 
dissociation probability contour plot. These are interpreted as remnants of the ellipticity of 
the bifurcated O e . We notice that the corresponding hyperbolic periodic orbit Oh remains 
hyperbolic (R < 0) for most values of the parameters. However, for 6 [tt/6, it], there is 
a region around the upper-right corner of the (Fi,F 2 ) plane where the Oh turns to elliptic 
and then returns to hyperbolic, as it is shown on the (Fi,F 2 ) plane contour plot of Fig. H] 
(6). In general, this bifurcation does not affect the dissociation probability because, due to 
its location, it does not play an important role compared to the periodic orbit O e which has 
already bifurcated (R > 1) for these parameter values (see Figs. H] (c) and (<f)). However, 
for low values of and high values of the amplitudes Fi, the orbit O e is still elliptic and 
Oh undergoes a bifurcation, as shown in Figs. [7] which we discuss in the following section. 
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This region corresponds to the disagreement observed for <fi = on Fig. H] with the direct 
simulations of Ref. 7j to which we turn next. 



C. Low values of eft : Influence of Oh 

For = 0, we observe two branches on Fig. H] (a) where the residues of O e vanish. Along 
these lines we expect, from a linear stability analysis, locally a constant degree of chaos, 
and hence dissociation probability. However, this is not seen in direct simulations since the 
dissociation probability actually increases along these lines as F 2 is increased. This feature 
is shown using laminar plots which represent contour plots of the number of return times 
on the Poincare section before dissociation (defined as trajectories for which E becomes 
greater than E t h = 2). The maximum integration time is 2007r/u;i ~ 2244. Figure [5] shows 
laminar plots for the set of parameters on the line of vanishing residues (as F 2 is increased) 
as marked by crosses in Fig. H] (a), and Figure [6] represents two additional laminar plots 
for parameter sets transverse to that line (with a fixed F 2 and increasing Fi), marked by 
circles in Fig. H] (a). These plots show clearly that dissociation increases as F 2 is increased. 
Even locally around the considered periodic orbit O e , there seems to be more chaos as F 2 
is increased. The situation is more complex at parameter values off the vanishing-residue 
line as Fx is varied on Fig. [H] where the overall amount of dissociation seems to be similar 
but with very different distributions of dissociating trajectories. All these features originate 
from the nonlinear stability which can be captured by considering the linear stability of 
higher order periodic orbits around the boundary of the elliptic island O e . 

For insights into the parameter region where this discrepancy occurs, we describe here the 
associated bifurcation. In Fig. [7] the residue as a function of F\ is plotted for both O e and 
Oh with fixed F 2 and 0. We see clearly that there is a loop around R = 0, indicating that 
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Oh undergoes a bifurcation at R = which involves three periodic orbits of the same period, 
two hyperbolic ones and an elliptic one. When is equal to 0, the upper part of the loop 
merges with the residue curve of O e , as shown in Fig. 0(6), for which the loop size reaches 
its maximum size. The following picture emerges : Without a loop in the residue curve, the 
system has two periodic orbits with the period of the field, O e and Oh- In the region of 
parameters where there is a loop in the residues, the system has four of these periodic orbits, 
two elliptic ones (which are close to each other or even coincide at <fi = 0) and two hyperbolic 
ones. This additional hyperbolicity increases chaos (and hence dissociation) locally. 

Conclusions 

We have analyzed the dissociation dynamics of a model diatomic molecule driven by 
a bichromatic field in terms of periodic orbit bifurcations. Following the linear stability 
of a few selected periodic orbits, we reproduced the dissociation probability qualitatively 
in parameter space (two field amplitudes and one relative phase). For relatively low 
and high amplitudes F i: the original hyperbolic periodic orbit Oh undergoes a particular 
bifurcation, which leads to two branch lines on the F1-F2 residue plane. Along these two 
branch lines, there is a discrepancy between predictions based on the residues and direct 
simulations. The role of additional periodic orbits is underlined regardless of whether the 
discrepancy originates from bifurcated orbits (and the resulted increase of hyperbolicity) or 
from higher-order periodic ones. 
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List of Figures 

1 Poincare section of Hamiltonian (J2J) for F\ = 0.18, F 2 = 0.02 and = 0. The 
dot and the cross indicate the elliptic periodic orbit O e and its associated 
hyperbolic one, respectively. The inset depicts a projection of O e in the (r,p) 
plane 

2 Action(solid) and angle(dashed) curves for F\ = 0.18 and F 2 = 0.02 

3 Residue curves (bold solid) for the considered set of periodic orbits for F\ = 
0.18 and F 2 = 0.02 : for O e (upper curve) and for Oh (lower curve). The 
arrow indicates where the bifurcation happens 

4 Contour plots of the residues of O e in the (Fi,F 2 ) plane for (a) = 0, 
(6) (p = 7r/6, (c) <fi = 7r/2 and (d) (p = ix. White dashed curves indicate 
R(Fi,F 2 ) = 1 where O e bifurcates. Two white solid curves in (a) indicate 
R(F\,F 2 ) = 0. The small region enclosed by two thin white lines near the 
up-right corner in (b) corresponds to a region where Oh is elliptic. Crosses 
(respectively circles) on panel (a) represent the values of parameters where 
laminar plots have been computed in Fig. \5\ (respectively Fig. E]) 

5 Laminar plots with <p = for (a) {F 1 ,F 2 ) = (0.1405,0.025), (b) (F u F 2 ) = 
(0.1559,0.04) and (c) (F 1 ,F 2 ) = (0.165,0.055). The cutoff time is 200^/^1 w 
2244 and the energy threshold is E t h = 2 

6 Laminar plots with = for (a) (F 1 ,F 2 ) = (0.19,0.055), and (6) (F 1 ,F 2 ) = 
(0.14,0.055). The cutoff time is 200^/^! w 2244 and the energy threshold is 
E t h = 2 
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7 Residue versus Fi for both elliptic and hyperbolic periodic orbits at F 2 = 0.03 
for (a) several different values of : <p = tt/6 (dash-dotted), = 0.35 (solid), 
= 0.2 (dotted) and = 0.06 (dashed), and (6) = [l9 



Figures 




FIG. 1: Poincare section of Hamiltonian ([2]) for Fx = 0.18, F 2 = 0.02 and = 0. The dot and the 
cross indicate the elliptic periodic orbit O e and its associated hyperbolic one, respectively. The 
inset depicts a projection of O e in the (r,p) plane. 
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15 




FIG. 3: Residue curves (bold solid) for the considered set of periodic orbits for F\ = 0.18 and 
F2 = 0.02 : for O e (upper curve) and for Oh (lower curve). The arrow indicates where the 
bifurcation happens. 
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FIG. 4: Contour plots of the residues of O e in the (Fi,^) plane for (a) <j> = 0, (6) <j> = n/6, (c) 
(j) = 7r/2 and (d) cp = tt. White dashed curves indicate R(Fi,F2) = 1 where O e bifurcates. Two 
white solid curves in (a) indicate R{F\ } F2) = 0. The small region enclosed by two thin white lines 
near the up-right corner in {b) corresponds to a region where Oh is elliptic. Crosses (respectively 
circles) on panel (a) represent the values of parameters where laminar plots have been computed 
in Fig. [5] (respectively Fig. [6]). 
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FIG. 5: Laminar plots with = for (a) {F 1 ,F 2 ) = (0.1405,0.025), (b) (F 1 ,F 2 ) = (0.1559,0.04) 
and (c) (F 1 ,F 2 ) = (0.165,0.055). The cutoff time is 200^/^1 « 2244 and the energy threshold is 
E t h = 2. 
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FIG. 6: Laminar plots with <f> = for (a) (-Fi,^) = (0.19,0.055), and (6) (F 1 ,F 2 ) = (0.14,0.055). 
The cutoff time is 200it/uji 2244 and the energy threshold is E t h = 2. 
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FIG. 7: Residue versus F\ for both elliptic and hyperbolic periodic orbits at F 2 = 0.03 for (a) 



several different values of 



7r/6 (dash-dotted), 



0.35 (solid), 



0.2 (dotted) and 



^ = 0.06 (dashed), and (6) <j) = 0. 
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